
# library(pacman)
# p_load(dplyr, ggplot2, stringr, tidyr, countrycode, DBI, RSQLite) 
library(dplyr)
library(tidyr)
library(stringr)
library(ggplot2)
library(countrycode)
library(DBI)
library(RSQLite)

options(scipen = 999) # Not scientific numbering on plot

con <- dbConnect(RSQLite::SQLite(), "iso_standards.sqlite") # Attaching to database

country_certifications <- dbReadTable(con, "country_certifications") # Extracting the certifications dataset

dbDisconnect(con) # Disconnecting from database

cbPalette <- c("#164688", "#f53f00", "#fad200", "#5f9d03", "#77001f", "#8ecbff",
               "#324000", "#b0cf00", "#492071")

country_certifications <- country_certifications %>% # Fixing some name errors (should be fixed in the current database version of the dataset)
  mutate(country = ifelse(country == "Russion Federation", "Russian Federation", 
                          ifelse(country == "Ethopia", "Ethiopia",
                                 ifelse(country == "Albana", "Albania",
                                        ifelse(country == "Taiwan, Province of China", "Taiwan",
                                               ifelse(country == "Taipei, Chinese", "Taipei",
                                                      ifelse(country == "Cocos (Keeling) Islands, Kosovo", "Kosovo",
                                                             ifelse(country == "Netherlands Antilles (NL)", "Netherlands Antilles",
                                                                    ifelse(country == "United Kingdom of Great Britain and Northern Ireland", "United Kingdom",
                                                                           country)))))))))


cert_agg <- country_certifications %>%
  mutate(continent = countrycode(sourcevar = country,
                                 origin = "country.name",
                                 destination = "un.region.name")) %>% # Standardizing names using countrycode
  mutate(iso_name = str_to_upper(iso), 
         iso_name = str_replace_all(iso_name, "_", " "),
         iso_name = str_squish(iso_name),
         iso_name = ifelse(iso_name == "ISO IEC 27001", "ISO 27001", 
                           ifelse(iso_name == "ISO IEC 20000-1", "ISO 20000 1", iso_name))) %>% # Getting all ISO series names into the same version
  mutate(certificates = replace_na(certificates, 0)) # NA indicates no certifications for that country in given year, replacing with 0

cert_agg %>%
  na.omit() %>%
  mutate(continent = ifelse(continent == "Asia", "Asia & Oceania",
                            ifelse(continent == "Oceania", "Asia & Oceania", continent))) %>% # Merging Asia and Oceania to one continent
  group_by(continent, year, iso, iso_name) %>% 
  summarise(certificates = sum(certificates, na.rm = TRUE)) %>% # Sum up number of certificates per continent year-wise
  filter(iso %in% c("iso_9001", "iso_14001", "iso_13485", "iso_iec_27001")) %>% # Filter out selected ISO series
  mutate(Standard = ifelse(iso == "iso_9001", "ISO 9001: Quality Management", # Making nice names out of the ISO series
                           ifelse(iso == "iso_14001", "ISO 14001: Environmental Management Systems",
                                  ifelse(iso == "iso_13485", "ISO 13485: Medical Devices",
                                         ifelse(iso == "iso_iec_27001", "ISO 27001: Information Security Management Systems"))))) %>%
  mutate(year = as.numeric(year)) %>% # Numeric for plot-friendly year axis
  mutate(` ` = continent) %>% # No need for label name on grouping variable
  ggplot(aes(year, certificates, group = ` `, color = ` `)) + # Plotting
  geom_line(aes(linetype = ` `), linewidth = 1) +
  scale_fill_manual(values=c(cbPalette), name = "") +
  facet_wrap(~ Standard, scales = "free") + 
  theme_classic() +
  labs(x = "", y = "Number of certificates") +
  theme(axis.title.y = element_text(size=rel(1.8)),
        legend.text=element_text(size=rel(1.8)),
        legend.position = "bottom",
        legend.title = element_text(size = rel(1.8)),
        axis.text.x = element_text(size = rel(1.8), angle = 90, vjust = 0.5, hjust=1),
        axis.text.y = element_text(size = rel(1.8)),
        strip.text = element_text(size = 14))
